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Abstract As large-distance rays (say, 10-24°) approach the solar surface ap- 
proximately vertically, travel times measured from surface pairs for these large 
separations are mostly sensitive to vertical flows, at least for shallow flows within 
a few Mm of the solar surface. All previous analyses of supergranulation have 
used smaller separations and have been hampered by the difficulty of separating 
the horizontal and vertical flow components. We find that the large separation 
travel times associated with supergranulation cannot be studied using the stan- 
dard phase-speed filters of time-distance helioseismology These filters, whose 
use is based upon a refractive model of the perturbations, reduce the resultant 
travel time signal by at least an order of magnitude at some distances. More 
effective filters are derived. Modeling suggests that the center-annulus travel 
time difference [6t i] in the separation range A = 10-24° is insensitive to the 
horizontally diverging flow from the centers of the supergranules and should lead 
to a constant signal from the vertical flow. Our measurement of this quantity, 
5.1 ± 0.1 seconds, is constant over the distance range. This magnitude of signal 
cannot be caused by the level of upfiow at cell centers seen at the photosphere of 
10 ms" 1 extended in depth. It requires the vertical flow to increase with depth. 
A simple Gaussian model of the increase with depth implies a peak upward flow 
of 240 ins -1 at a depth of 2.3 Mm and a peak horizontal flow of 700 ms -1 at a 
depth of 1.6 Mm. 
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1. Introduction 



Supcrgranulation, first studied more than 50 years ago QHart, 19 54 Lcig hton, Noyes, a nd Simon, 1962), 



continues to be an active area of research. A comprehensive review ( |Rieutord and Rincon, 2010[ )[ 
details the recent progress. It has proven difficult to measure the vertical flow of 
supergranules at the photospheric level. The recent measurements of a 10 ms -1 
upflow at the center of the average cell with a horizontal variation consis- 
tent with a simple convection cell flDuvall and Birch, 2010[ ) may have finally 
settled the issue. It is only through helioseismology that we would hope to mea- 
sure the supergranulation flows below the photosphere. IDuvall et al.\ (|1997[) first 
showed that time-distance helioseismic techniques are sensitive to supergranules 
and that inversions to derive three-dimensional flow fields might be derived. 
However, I Zhao and Kosovichevl (|2003p showed that their ray-theory inversions 
could not separate the horizontal and vertical flows for a model flow field. 
The use of radial-order filters and Born approximation kernels has led to more 



successful separation flJackiewicz, Gizon, and Birch, 2008 Svanda et at, 2011 1 



|Gizon, Birch, and Spruft] (|2010p summarize the local helioseismic contributions 
to supergranulation. 

An important quantity in time-distance helioseisomology is the arc separa- 
tion [A] between pairs of photospheric locations whose signals are subsequently 
temporally cross correlated. Previous analyses used A less than 5° (61 Mm) 
([Zhao and Kosovichev,~2003| or 2.4° (29.2 Mm) flJackiewicz, Gizon, and Birch, 2008) 1 



Svanda et al, 2011 ), as the sensitivity is greater for the horizontal flow, the 
signal-to-noise ratio is better for small separations, and basically no super- 
granulation signal could be seen at larger A in the travel time difference maps 
constructed from the relatively short (8-12 hours) time intervals required to 
study the one-day lifetime supergranules. In the present work, we show that the 
large separations of A = 10-24° yield center-annulus travel time differences 
[St Q i = outward-going time minus inward going time] from the centers of average 
supergranules that are insensitive to the diverging horizontal flow and hence yield 
a purely vertical flow response. 

In Section [2] the modeling efforts, including flow models that satisfy mass 
conservation, travel times computed from ray theory, and linear wave simulations 
through flow fields are presented. In Section [21 the analysis of data is presented, 
and in Section @] the results are discussed. 



2. Modeling 
2.1. Flow Model 

As our observational technique is to make averages about centers of super- 
granules, the flow model that we take is azimuthally symmetric and decays 
exponentially away from the cell center to simulate the effect of averaging a 
stochastic flow field. A time-independent flow is assumed with axisymmctry 
about the location of the average cell. A Cartesian coordinate system with z 
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positive upwards is denned. For these assumptions, the vanishing divergence of 
the mass flux implies 

d 

Q^(P V *) = -P^h ■ v h , (1) 

where p is the z-dependent density, v z is the z-component of velocity, Vh is 
the horizontal divergence, and Vh is the horizontal velocity. It is assumed that 
horizontal density variations can be ignored, as travel times are only weakly 



sensitive to density changes in the underlying medium (Hanasoge et al, 20121 



The model is assumed to separate into horizontal and vertical functions 

v h = -f(z) S (x,y), (2) 



v z = u(z)V h • g, (3) 

where g(x,y) is a vector in the horizontal plane with no units, /(z) has units 
of velocity, and u{z) has units of velocity times length. Our method of solution 
is to first specify the horizontal function g(x, y) and calculate its horizontal 
divergence Vh • g- The vertical function u{z) is then specified, and Equation (p} 
is used to derive f(z). Some straightforward algebra yields 

n / \ du din/? 

fi*) = Tz +u W-aT- (4) 

In general, models are considered with the horizontal function g(x, y) defined by 

g (x,y)=rJ 1 (kr)e- r / R , (5) 

where f is the outward radial unit vector in a cylindrical coordinate system, J\ is 
the Bessel function of order one, fc is a wavenumber, r is the horizontal distance 
from the origin, and R is a decay length, k and R are, in principle, free parameters 
that will only have the default values k = |^ rad Mm -1 and R = 15 Mm in this 
article. This type of horizontal variation was used by IBirch and Gizonl (|2007p . 
This type of horizontal variation describes well the type of averaging of the 
supcrgranular field done in this paper, where one sees the outward flow as the 
first positive lobe of the J\ function and on average, the inflow from the adjacent 
supergranular cells as the first negative range of the J\ function. 

In general, models have been considered with a Gaussian z-dependence. v z is 
specified at r = as a simple Gaussian: 

v z {r = 0) = ku{z) = v e~ (z - zo)2/2c7 ', (6) 

where zq is the location of the peak of the vertical flow, a z is the Gaussian sigma, 
and vq is the maximum vertical flow. To ultimately explain the large-distance 
travel times, the photosphere needs to be in the far tail of the Gaussian. As 
the upward flow at cell center is 10 ms -1 at the photosphere, this implies a 
considerably larger vertical flow at depth for the average cell. Values for these 
parameters that approximate the data are vq = 240 ms -1 , zq = —2.3 Mm, and 
a z = 0.912 Mm. 
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Figure 1. Comparison of various ray models. The blue curves are for the Gaussian model 
with zo = —2.3 Mm. The lower dot-dashed blue curve is for the horizontal flow component, the 
upper dot-dashed blue curve is for the vertical flow component and the solid blue curve is for 
the sum (upper solid curve). The upper dashed (green) curve is for the sum for zo = —1.15 Mm 
model and the lower dashed (red) curve is for the sum of the zq = —3.45 Mm model. The (lower) 
solid black curve is the sum for a model with a constant upflow of 10ms -1 . 



2.2. Ray Calculations 

In the ray theory the travel time difference for the two directions of propagation 
through a flow v is 

/ r ^, (7) 

where T is the ray path, c is the sound speed, and ds is the element of length in 
the direction of the ray ( |Giles, 2000[ ) . This equation naturally separates into two 
terms for the horizontal and vertical flow contributions to the travel time. The 
ray generation and raw ray integrations are performed using code developed 
and discussed in detail by IGilesI (|2000|) . This code was extended to integrate 
quadrant and annulus surface-focusing integrations for 3D flow and sound speed 
perturbation models. In Figure[T]the black curve shows the travel time difference 
[St i] for a flow model with a constant upflow of 10ms" 1 . The black curve is 
for the sum of horizontal and vertical flows. That is, a horizontal flow of the 
form of Equation ((5|) consistent through the continuity equation with a 10 ms -1 
upflow is used. This signal, with a maximum of 0.6 seconds, is much too small 
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Table 1. Model details. The three models gl, g2, g3 are ones with Gaussian (g) 
upward flow at cell center peaking at zOMm with peak flow v z (zO) ms~ 1 and with 
Gaussian width cr z Mm. The upward flow at the top (z=0) is v z (0)ms~ 1 . The 
peak horizontal flow at the top is ti™ ax (0)ms which occurs at radial distance 
r max Mm. The peak outward flow «™ ax ms _1 occurs at 2(u™ ax )Mm. All of the 
maxima of Vh occur at the same r max by construction. The model siml is the one 
for the small features with just a z flow that are used for the flow perturbation 
described in Section 12.31 



Name 


zo 




vz(0) 


» 4 mal (0) 


j,max 


v z (zO) 


..max 
v h 


2(>™ aX ) 


gl 


-3.45 


1.39 


10 


122 


7 


218 


460 


-2.44 


g2 


-2.3 


0.912 


10 


138 


7 


240 


697 


-1.62 


g3 


-1.15 


0.433 


10 


195 


7 


340 


1609 


-0.83 


siml 


-2.3 


0.816 


6.3 






338 







to explain the results derived in Section 13.31 This implies that the vertical flow 
must increase with depth from the photospheric value of 10 mr 1 . 

Also in Figure [TJ the separate horizontal and vertical travel time contribu- 
tions and the sum are shown for the Gaussian model mentioned above. For this 
relatively shallow model, it is seen that in the distance range 10-24°, the signal 
is essentially constant and is mostly due to the vertical flow. Towards shorter 
distances, the horizontal flow contributes an increasing amount. The Gaussian 
model has three free parameters: vq, zq, and a z . In a later section the travel time 
difference for 10-24° will be determined to be 5.1s. This travel time difference 
and the 10 ms -1 upflow at the photosphere serve to determine two of the three 
free parameters. The travel time at A = 24.5° is forced to 5.1 seconds. Letting 
zq be the free parameter, three models are plotted in Figure [JJ with values of 
z = -1.15 Mm, z = -2.30 Mm, and z a = -3.45 Mm. The shallowest model 
shows the least contribution from the horizontal component in the distance range 
5- 10°. This distance range may supply a way to distinguish the best model while 
comparing with data. 

A summary of the characteristics of these models is contained in Table [TJ The 
maximum vertical flow is at least a factor of 20 larger than the photospheric 
value. The ratio of the maximum horizontal flow to the photospheric value varies 
between four and eight. The maxima of the horizontal flow occurs somewhat 
shallower than the peak vertical flow. 

2.3. Simulation Results 

The applicability of ray theory to time-distance helioseismology has been called 

into question ( |Bogdan, 1997D . There are cases for supergranulation studies (jBirch and Gizonl | 

120071 Figure 5) where it works reasonably well (amplitude of ray theory 25% 

higher than Born approximation kernels) and others where it is discrepant 

by a factor of two (jBirch and Gizonl 120071 Figure 6). One case that works 

extremely well is the comparison of interior rotation determined from global 

modes and time-distance inversions ([Giles , 2000. Figure 6.2) in the radius range 

0.89 < r/ii© < 0.999. The difference between these cases seems to be that the ray 



SOLA: ms.tex; 26 July 2012; 0:29; p. 5 



T.L.Duvall Jr.,S.M.Hanasoge 



8 



S 4 




A [deg] 

Figure 2. Comparison of the center-annulus travel time differences [<5t ;] from the linear 
simulation (blue) with the travel-time difference computed from the ray theory with the 
same flow perturbations (red). The error bars are computed from the scatter far from the 
feature locations. No filtering has been done before the travel time measurements. The average 
travel-time difference in the range A = 10-24° has been scaled to match the observationally 
determined mean 5.1 seconds. The same scaling factor is then used to scale the ray theory 
results. In the range A = 10 — 24°, the ray theory predicts a travel— time difference too large 
by 24%. 



theory works better when the spatial variation of the perturbation is larger, that 
is, much larger than the acoustic wavelength. One way to check the ray theory 
for a particular case is to compare with travel times computed from Born kernels 
(Bir ch and Gizon, 2007[ ). Another way to check is to perform simulations on a 
convectively-stabilized solar model with specified flow velocity perturbations and 
to propagate acoustic waves through the model ( |Hanasoge, Duvall, and Couvidat, 2007| ) .| 
Travel times are computed from the simulations and then compared with ray- 
theory calculations of the travel times using the same flow velocity perturbations. 
This is the type of checking of the ray theory adopted here. 

A global simulation of wave propagation ( Hanasoge and Duvall, 2007 1 with 
768 x 384 grid points in longitude and latitude respectively is performed over 
a ten-hour period. The flow perturbation consists of Gaussian features with 
a constant direction of vertical flow centered at a depth of zq = —2.30 Mm 
with a z = 0.82 Mm and with horizontal er^ = 5.1 Mm. There are 500 of these 
features placed at random longitude-latitude pairs at that depth [zq]. Center- 
annulus travel time differences arc measured from the simulation and averaged 
about the known locations. Realization noise is removed to first order by doing a 
similar simulation with no flow perturbation but with the same source excitations 
QWerne, Birch, and Julien,~2 004; Hanaso ge, Duvall, and Couvidat, 2007| ) and sub-| 
tracting the resultant 5t j to obtain noise-corrected results. 

The flow model that is inserted into the simulation is used with Equation ([7]) 
to derive ray-theory predictions of the travel time differences. The results of 
the ray-theory computations and the travel time difference measurements from 
the simulation are shown in Figure [5] As the simulation only contains vertical 
flows, it is expected that there would be little or no variation of the travel-time 
differences with A. The ray theory does predict too large a travel time difference 
by 24%. This excess amplitude is similar to that found by IBirch and Gizonl 
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(2007). Although the amplitude computed from the ray theory is not precise, the 
necessity of large vertical flows to generate a travel-time difference of 5.1 seconds 
at A = 10-24° is confirmed. The peak vertical flow for the normalized case is 
338 ms -1 . The model parameters are detailed in Table [TJ 



3. Data Reduction and Analysis 

3.1. Reduction 

Dopplcrgrams from the Helioseismic and Magnetic Imager (HMI) : ( |Schou et al., 2012[ )[ 
onboard the Solar Dynamics Observatiory (SDO) spacecraft were analyzed for 
the present work. 32 days (10 June- 11 July 2010) of Dopplcrgrams were used 
to derive the final results, although a number of tests subsequently mentioned 
used only the last three days, 9-11 July 2010. This particular time period was 
used as the Sun was very quiet (sunspot number RI w 15 for these two months), 
and it was the final two-month continuous coverage period for the Michelson 
Doppler Imager (MDI) instrument (Scherre r et al., 1995| ). It might be useful to 



compare the results from the two instruments, but only the HMI data was used 
for the present study. 

Raw Dopplergrams have a nine-day averaged Dopplergram subtracted as well 
as the spacecraft velocity which has a significant 24-hour component. These 
corrected Dopplergrams are remapped onto a coordinate system with equal 
spacing in latitude and longitude of 0.03° covering a range of 144° in both 
latitude and longitude. The remapping is achieved using a bilinear interpolation. 
The remapping in longitude is onto the Carrington system at a central longitude 
that crosses central meridian at the middle of the twelve-hour interval used. Two 
twelve-hour intervals are used for each day, covering the first and second halves 
of the day. Each remapped image is Fourier-filtered and resampled at 0.24° 
per pixel. Twelve-hour datacubes are constructed from these individual images 
with 600 x 600 spatial points and 960 temporal points for the 45 s sampling. 
This procedure works well for the large A emphasized in the present study, but 
would need to be modified to study smaller A. Studying smaller A is deferred 
until future work. 

The centers of the supergranules are located by a procedure used previously 



(Du vall and Birch, 2010 ). The datacubes are spatio-temporally filtered to pass 
just the solar /-mode oscillations. Center-annulus travel time differences are 
computed ( |Duvall and Gizon, 2000[ ) for the distance range 0.48-1.02°. By us- 
ing the difference of inward-going times minus outward-going times, the maps 
are equivalent to maps of the horizontal divergence. The equivalence of these 
maps to maps of supergranules has been shown before ( |Duvall and Gizon7"2 000; 



Duvall and Birch, 2010 1. The travel-time difference maps arc smoothed with a 
Gaussian of a = 2.9 Mm to more easily determine the cell-center locations. Local 
maxima of this signal are picked as the cell centers. Lists of the locations are 
made. The one of a pair with the smaller divergence that are closer together 
23 Mm is rejected. As it is desired to use center-annulus separations A <= 25°, 
locations close to the edge of the maps are not used. On average, 930 features 



SOLA: ms.tex; 26 July 2012; 0:29; p. 7 



T.L.Duvall Jr.,S.M.Hanasoge 



Table 2. Nominal phase speed filter parameters. The first col- 
umn is an identification number. The second column is the range 
in distance [A]. The third column is the mean distance [A mea n]- 
The fourth column is the central phase speed. The fifth column is 
the full width at half maximum (FWHM) of the Gaussian phase 
speed filter. For the tests in Section 13.21 filters 1-13 were used 
while for the final results in Section 13.31 filters 1-14 were used. 



Filter 


A[°] 


A mean [ ] 


^0 [ km s 1 ] 


<5u[kms - 1 ] 


1 


1.56-2.28 


1.92 


32.0 


13.2 


2 


2.28-3.00 


2.64 


41.8 


4.2 


3 


3.00-3.96 


3.48 


46.9 


6.0 


■1 


3.96-4.92 


4.44 


54.7 


9.5 


5 


4.92-6.12 


5.52 


64.5 


10.1 


6 


6.12-7.80 


6.96 


75.1 


11.1 


7 


7.56-9.24 


8.40 


84.3 


10.4 


8 


9.12-10.7 


9.84 


93.5 


9.4 


9 


10.92-12.60 


11.76 


104.4 


9.7 


10 


12.36-14.52 


13.44 


114.2 


12.5 


11 


14.28-16.44 


15.36 


125.2 


12.3 


12 


17.16-19.32 


18.24 


141.9 


12.6 


13 


19.08-22.20 


20.64 


156.0 


18.5 


14 


22.08-24.00 


23.04 


170.4 


11.8 



are located for each datacube with a total of 59549 for the 64 twelve-hour 
datasets. Overlays of cell-center locations with the divergence maps are shown 



in (Duvall and Birch, 2010). 



3.2. Filtering 

Ray theory has proven invaluable in the development of time-distance helioseis- 
mology ( |Duvall et ai, 1993] ) . It led to the basic idea that signals propagating 
along ray paths between surface locations would lead to correlations between 
the temporal signals from which travel times could be inferred. As acoustic, or 
p, modes with the same horizontal phase speed \uj/k: tu angular frequency and 
k horizontal wavenumber] travel along the same ray path to first order, it was 
natural to filter the data in ui/k to isolate waves traveling to a certain depth. This 
type of filter was introduced by IDuvall et al.\ jTT>').''»! with subsequent develop- 
ment bv lKosovichev and Duvalll(|1997p and [Giles ( 20DS1)- I n the spatio-temporal 
power spectrum of solar oscillations {k~uj diagram), a line from the origin has a 
constant phase speed uj/k. A range of phase speeds thus corresponds to a pie- 
shaped wedge. Phase-speed filtering then corrcponds to multiplying the Fourier 
transform of the input datacube by a pie-shaped wedge (appropriately modified 



for 3D) with some Gaussian tapering in phase speed ( |Couvidat and Birch, 2006 1 . 

The overriding principle of the use of the phase-speed filter is that solar per- 
turbations lead to refractive time changes for signals travelling along a ray path. 
That this model is deficient was shown nicely in the work offCouvidat and Birchl 
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Figure 3. Center-annulus travel time differences [<54 ;] averaged about the supergranulc 
centers for different filters for the thirteen A ranges for three days of HM1 data. From the 
bottom to top (on the right side), the black (x) points are for the the nominal phase-speed 
filters. The blue (*) points are for the phase-speed filters with the widths doubled. The red (o) 
points arc for the width tripled. The orange (square) points are for a filter width that is half the 
nominal phase speed. The magenta (diamonds) curve is for the constant dcgrec-width-filtcr 
(Section 13. 311 with width Vi = 400. The green (pentagrams) curve is for no phase-speed filter. 
In all cases the / mode is excluded via filtering as is signal outside the frequency bandpass 
1.5 < v < 6mHz. These results were obtained using the three days of data 9-11 July 2010. 

(2006) in which an artificial signal was detected with different amplitude de- 
pending on the width [5v] of the pie-shaped wedge. They concluded that the 
Born-approximation kernels ( |Birch and Kosovichev, 2000[ ), which are derived 
based on a single scattering from solar perturbations, arc more appropriate. 
A better way to think of the problem is to consider waves impinging on a 
perturbation which are subsequently scattered. 

An initial analysis was performed using some nominal phase speed filters and 
distance ranges detailed in Table [2J Thirteen distance ranges (and hence thirteen 
different filters) were applied to the six input datacubes for the dates 9-11 July 
2010. The overall distance range covered for this test is 3.0-22.2°. Center- 
annulus travel-time differences [<5i i] were obtained for the thirteen distance 
ranges for each of the six datacubes using a standard surface-focusing time- 
distance analysis with the Gizon-Birch method of extracting travel times from 
the correlations ( |Gizon and Birch, 2004[ ). The travel times are averaged over 
the distance range and about the locations of the supcrgranular centers and the 
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signal at the average cell center location is extracted. Uncertainties are estimated 
from the scatter of the six values. The travel-time differences are plotted versus 
A in Figure [3] (black). The variation with A does not agree with the models in 
Figurc[TJ but peaks near A = 5° and decays approximately to zero near A = 20°. 

There was some concern that this distance dependence might have something 
to do with the filtering, so a case was analyzed with no phase-speed filtering. 
The / mode was excluded and there was a frequency bandpass filter transmitting 
v = 1.5-6mHz. The results of that analysis are also shown in Figure [3] (green). 
This shows a behavior with A much closer to the models in Figure [1] but with 
considerably more noise and a much larger signal close to five seconds. It was 
concluded that the phase speed filters were significantly degrading the signal. 
To test this further, the filters were broadened by the factors 2, 3, 4, 5 with the 
results of factors 2 (blue) and 3 (red) also shown in Figure [3] More signal is 
obtained with the larger-width filters, but even with the factor of 5, only about 
50% of the signal is obtained and there is still a decay of the signal towards 
larger A. 

Because of the approximate constancy of the simulation travel times [St a i] 
with A in Section 12.31 and the approximate constancy with A for the no phase- 
speed filter case, it was concluded that the decay of 5t Q i toward larger A in 
the phase-speed filtered results is an artifact of the filtering. Filtering would 
seem to be desirable because of the large difference in errors for the unfiltered 
and filtered cases. But a filter is needed that gives a more unbiased result at the 
large A in order that the results not be overly dependent on the modeling. A first 
attempt at this was to use a set of filters with the full width at half maximum 
(FWHM) equal to half the phase speed ( |Birch et ai, 2006 ) . The results of these 
measurements are also shown in Figure [3] (orange). The dependence on A is 
somewhat reduced but there is still some decay at the largest A. 

A different way to look at this problem is to consider a p mode with frequency 
v and spherical harmonic degree £ impinging on a supcrgranulc. To the extent 
that the supcrgranulc is static, the mode frequency v will be maintained and 
the scattering will spread power in the spectrum away from the nominal value 
of I. The scattering should spread power over a width in I that depends on 
the spectrum of supcrgranulation, which peaks near I = 120. For similar values 
of the phase speed (which can also be characterized by v/i), the spread power 
should be over a similar width in i QChou et al., 1996} |Woodard, 2002[ ). A filter 
that is centered at a certain phase speed but whose width is constant with v 
and characterized by the FWHM in I (= Ti) may be what is needed. A value 
of Yi — 400 should capture most of the supergranular signal. Such a filter was 
implemented and applied to the thirteen distance ranges used. The results are 
also shown in Figure [3J The approximate constancy with A > 5° and the larger 
values than for the normal phase speed filters suggests that this type of filtering 
is a useful concept. The filter used has a flat top of width T^/2 and is tapered 
to zero by a cosine bell that goes from 1 to in Tg/2. 

Some additional tests of this new filter concept were computed. In Figure[UJa) 
is shown the peak supergranular signal for A = 19.08-22.2° as a function of the 
filter width Ti. For small Ti, much of the supergranular signal is removed, but 
as the width is increased more and more of the signal is transmitted. There 
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Figure 4. (a) Travel-time difference [<5t i] versus filter FWHM [I^]. The unfiltered case has 
5.3 ± 1.2 seconds, (b) The travel time difference from (a) divided by the size of the error bar 
from (a) versus the filter FWHM Tp. The value for the unfiltered case is 4.6. 



is clearly a variation in the size of the error bars across Figure IHa). Another 
thing to examine is the ratio of the signal in Figure EJa) to the size of the error 
bar, which should give a value of the signal to noise ratio (S/N) of the filter 
( Couvidat and Birch, 2006 ). This is plotted in Figure E£b). A rather broad peak 
is seen more or less centered on = 400. One noticeable aspect of Figure HJb) 
is that the unfiltered case gives a worse S/N (4.6) than the filters near the peak. 
Also narrow filters have reduced S/N. Based on the flatness of the Ti = 400 
results with A in Figure [3J the extraction of essentially all the supergranular 
signal in Figure Eta), and the value of the S/N in Figure E^b), Ti — 400 has 
been chosen for further analysis. 



3.3. T e = 400 Analysis 



The 64 twelve-hour datacubes were analyzed with the fourteen phase-speed 
filters with = 400. Center-annulus cross correlations were computed for 
the various distance ranges of Table 1. For checking purposes, it would be 
useful to be able to compute travel times using both the Gabor wavelet method 
dDuvall et al, 1997| ) and the Gizon-Birch method ( jGizon and Birch, 2004| |. How-| 
ever, for the large distances in the present study (A < 24°), individual twelve- 
hour cross correlations have insufficient signal-to-noise ratio to enable the Gabor 
wavelet fitting. It was decided to average the temporal correlations spatially 
about the supergranular centers. The averaging of the 930 (on average) corre- 
lations yields a high signal-to-noise ratio correlation that is amenable to the 
Gabor wavelet fitting at the largest distances [A]. Averaging the resultant <5i i 
for the 64 twelve-hour datacubes enables the use of the envelope-time differences 
[tenv] as well as the normal phase times [t p h]- Although most theoretical work 
applies to phase time differences, it will be useful to examine the envelope-time 
differences as they may have independent information. 

It was found previously that there is an offset from zero (of unknown origin) 
of the center-annulus travel-time differences ( |Duvall and Gizon, 2000[ ) . For the 
small A of that study, the value of the offset was 0.16 ± 0.02 seconds. This offset 
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Figure 5. (a) Center-annulus travel-time difference averaged about supergranule centers for 
the range of A = 22.08-24°. The scale of the colorbar at right is in seconds. Only the central 
point corresponding to the supergranular center is used in the present study, (b) Azimuthal 
average of (a). Note the offset at large radii. This offset is believed to be an artifact which 
needs to be removed from the results, (c) The offset at r = 58 Mm for the different travel-time 
definitions versus A. Blue is for the Gabor wavelet phase time differences. Red is for the 
Gizon-Birch phase time differences and the green is for the Gabor wavelet envelope time 
differences, (d) The resultant travel time differences averaged for the 64 12-hour datacubcs 
corrected for the offset in (c). The colors are the same as in (c). 



is assumed unrelated to the flows to be measured and needs to be removed 
from computed travel times. To measure it for the present study, the average 
maps about the supergranule centers for the 64 12-hour datacubes are com- 
puted for the different distance ranges. An example is shown for the largest 
range of A = 22.08-24° in Figure [UJa). An azimuthal average of Figure [5ja) 
is computed and the result is shown in Figure [SJb). At large radii [r], the the 
signal is approximately constant. This constant, which is assumed to be the 
offset that needs to be subtracted from the results, is plotted for the various A 
and travel-time methods in Figure 03c) . The Gizon-Birch and Gabor wavelet 
phase times yield almost exactly the same results, which is expected for the 
quiet Sun and the use of the same time window for fitting (20 minutes). The 
travel-time differences corrected for this offset are plotted in Figure E^d). Again, 
the Gizon-Birch and Gabor wavelet phase times yield almost precisely the same 
results with an average of 5.1 ±0.1 seconds on the range A = 10 - 24°. The Gabor 
wavelet envelope times yield a considerbly larger value of 9.7 ± 0.3 seconds on 
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Figure 6. Comparison of the three ray models from Figure^ with the HMI results (black 
with symbols and errors) of the = 400 filtering and the Gabor-wavelet phase— speed time 
differences from Figure (3d). These models are the sum of vertical and horizontal signals for 
the three Gaussian vertical flows peaking at zo = —3.45 Mm (red dot— dashed), zo = —2.3 Mm 
(blue solid), and zq = —1.15 Mm (green dashed) specified in Table [T] 



A = 10-24°. In general the group velocity, duj/dk is about one half of the 
phase speed u>/k for p modes, possibly leading to the larger travel time. The 
systematic error from Figure [5jc) is generally less the 10% of the signal shown 
in Figure (5Jd). 

The Gabor-wavelet phase-time differences are plotted in Figure [S] versus the 
three ray models of Figured] The increase of the observed times for A = 2.5-8° 
may be the effect of the horizontal flow. The agreement is a little better for the 
zo = —2.3 Mm model. 



4. Discussion 

The main observational result of this article is that in the distance range A = 10 - 
24° that the mean travel-time difference [5t Q i] at the center of the average 
supergranule is 5.1 ± 0.1 seconds. How secure is this result? There are no other 
results to compare with in this distance range. Also, the method of averaging over 
many supergranule s is uncommon ( |Birch et ai, 2006} |Duvall and Birch, 2010| 



Svanda et aL, 2011). The largest distance [A] used previously is 5° bv lZhao and Kosovichev[ 
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(|2003p . The phase-speed filter used in that study would have reduced the travel- 
time difference [<5t i] by a factor of five and so it is difficult to compare. In 
addition, only inversions were presented and not raw travel-time differences. 

One of the most interesting aspects of the present results is the large factor 
(> 20) between the photosphcric vertical flow and the peak vertical flow. This 
increase of vertical flow then also requires an increase of the horizontal flow from 
the photosphere to the peak of a factor between 3.7 and 8.3 from the simple 
Gaussian models gl-g3. Some simple tests suggest that this horizontal velocity 
should be detectible by a quadrant time-distance analysis, which will be left for 
future work. Some previous analyses have detected some larger flows below the 
surface than at the surface ( |Duvall et ai, 1997 ). Possibly the subsurface flow 
is highly variable which has made average properties difficult to ascertain from 
inversions of individual realizations. 

These results confirm that supergranulation is a shallow phenomenon, with 
the flow peaking within a few Mm of the photosphere. This would tend to 
support models in which supergranulation is related to a near-surface phe- 
nomenon such as granulation (Rast, 2003). The supergranular wave measure- 



ments ( Gizon, Duvall, and Schou, 2003 Sc 



rou, 2003 ), and in particular the vari- 



ation with latitude of the anisotropy of wave amplitude, would suggest a con- 
nection with rotation. The near-surface shear layer has been modeled as a 
possible exciter of supergranulation ( Green and Kosovichev, 2006[ ). |Hathaway| 



(1982) found that supergranulation as a convective phenomenon generates a 
near-surface shear layer. A difficulty with the near-surface shear layer as su- 
pergranular exciter is that recent simulations of IStein et al.\ (|2009[) develop a 
near-surface shear layer but no hint of excess power at supergranular scales is 
seen. 

Acknowledgements The data used here are courtesy of NASA/SDO and the HMI Science 
Team. We thank the HMI team members for their hard work. This work is supported by NASA 
SDO and the NASA SDO Science Center program through grant SCEX22011D awarded to 
NASA GSFC. S. M. H. acknowledges funding from NASA grant NNX11AB63G. 



References 

Birch, A., Duvall, T.L. Jr., Gizon, L., Jackiewicz, J.: 2006, Bull, of the Am. Astronom. Soc. 
38, 224. 

Birch, A.C., Gizon, L.: 2007, Astronomische Nachrichten 328, 228. 
Birch, A.C., Kosovichev, A.G.: 2000, Solar Phys. 192, 193. 
Bogdan, T.J.: 1997, Astrophys. J. 477, 475. 

Chou, D.-Y., Chou, H.-Y., Hsieh, Y.-C, Chen, C.-K.: 1996, Astrophys. J. 459, 792. 
Couvidat, S., Birch, A.C.: 2006, Solar Phys. 237, 229. 

Duvall, T.L. Jr., Jefferies, S.M., Harvey, J.W., Pomerantz, M.A.: 1993, Nature 362, 430. 
Duvall, T.L. Jr., Kosovichev, A.G., Scherrer, P.H., Bogart, R.S., Bush, R.I., deForest, C, 

Hoeksema, J.T., Schou, J., Saba, J.L.R., Tarbell, T.D., Title, A.M., Wolfson, C.J., Milford, 

P.N.: 1997, Solar Phys. 170, 63. 
Duvall, T.L. Jr., Birch, A.C.: 2010, Astrophys. J. Lett. 725, L47. 
Duvall, T.L. Jr., Gizon, L.: 2000, Solar Phys. 192, 177. 
Giles, P.M.: 2000, Ph.D. Thesis, Stanford Univ. 
Gizon, L., Duvall, T.L. Jr., Schou, J.: 2003, Nature 421, 43. 
Gizon, L., Birch, A.C.: 2004, Astrophys. J. 614, 472. 



SOLA: ms.tex; 26 July 2012; 0:29; p. 14 



Vertical Flows 



Gizon, L., Birch, A.C., Spruit, H.C.: 2010, Ann. Rev. of Astron. Astrophys. 48, 289. 

Green, C.A., Kosovichev, A.G.: 2006, Astrophys. J. 641, L77. 

Hanasoge, S.M., Duvall, T.L. Jr.: 2007, Astronomische Nachrichten 328, 319. 

Hanasoge, S.M., Duvall, T.L. Jr., Couvidat, S.: 2007, Astrophys. J. 664, 1234. 

Hanasoge, S.M., Birch, A.C., Gizon, L., Tromp, J.: 2012, Phys. Rev. Letters (submitted). 

Hart, A.B.: 1954, Mon. Not. Roy. Astron. Soc. 114, 17. 

Hathaway, D.H.: 1982, Solar Phys. 77, 341. 

Jackiewicz, J., Gizon, L., Birch, A.C.: 2008, Solar Phys. 251, 381. 

Kosovichev, A.G., Duvall, T.L. Jr.: 1997, editors F.P. Pijpcrs, J. Christcnscn-Dalsgaard, and 
C.S. Rosenthal, SCORe'96 : Solar Convection and Oscillations and their Relationship 
225, Kluwer, Dordrecht, 241. 

Leighton, R.B., Noyes, R.W., Simon, G.W.: 1962, Astrophys. J. 135, 474. 

Rast, M.P.: 2003, Astrophys. J. 597, 1200. 

Rieutord, M., Rincon, F.: 2010, Liv. Rev. Sol. Phys. 7, 2. 

Scherrer, P.H., Bogart, R.S., Bush, R.I., Hoeksema, J.T., Kosovichev, A.G., Schou, J., Rosen- 
berg, W., Springer, L., Tarbell, T.D., Title, A., Wolfson, C.J., Zayer, I., MDI Engineering 
Team: 1995, Solar Phys. 162, 129. 

Schou, J.: 2003, Astrophys. J. Lett. 596, L259. 

Schou, J., Scherrer, P.H., Bush, R.I., Wachter, R., Couvidat, S., Rabello-Soares, M.C., Bogart, 
R.S., Hoeksema, J.T., Liu, Y., Duvall, T.L. Jr., Akin, D.J., Allard, B.A., Miles, J.W., 
Rairden, R., Shine, R.A., Tarbell, T.D., Title, A.M., Wolfson, C.J., Elmore, D.F., Norton, 
A.A., Tomczyk, S.: 2012, Solar Phys. 275, 229. 

Stein, R.F., Georgobiani, D., Schafenberger, W., Nordlund, A., Benson, D.: 2009, AIP CS 
1094, 764. 

Svanda, M., Gizon, L., Hanasoge, S.M., Ustyugov, S.D.: 2011, Astron. Astroph. 530, A148. 
Werne, J., Birch, A., Julien, K.: 2004, editor: D. Danesy, In: SOHO 14 Helio- and 

Asteroseismology: Towards a Golden Future 559, ESA, Noordwijk, 172. 
Woodard, M.F.: 2002, Astrophys. J. 565, 634. 

Zhao, J., Kosovichev, A.G.: 2003, Editor: H. Sawaya-Lacoste, In: GONG+ 2002. Local and 
Global Helioseismology: the Present and Future ESA SP-517, ESA, Noordwijk, 417. 



SOLA: ms.tex; 26 July 2012; 0:29; p. 15 



SOLA: ms.tex; 26 July 2012; 0:29; 



Solar Physics 
DOI: 10.1007/* 



Subsurface Supergranular Vertical Flows as Measured 
Using Large Distance Separations in Time— Distance 
Helioseismology 

T.L. Duvall Jr. 1 • S.M. Hanasoge 2 ' 3 

© Springer •••• 

Abstract As large-distance rays (say, 10-24°) approach the solar surface ap- 
proximately vertically, travel times measured from surface pairs for these large 
separations are mostly sensitive to vertical flows, at least for shallow flows within 
a few Mm of the solar surface. All previous analyses of supergranulation have 
used smaller separations and have been hampered by the difficulty of separating 
the horizontal and vertical flow components. We find that the large separation 
travel times associated with supergranulation cannot be studied using the stan- 
dard phase-speed filters of time-distance helioseismology These filters, whose 
use is based upon a refractive model of the perturbations, reduce the resultant 
travel time signal by at least an order of magnitude at some distances. More 
effective filters are derived. Modeling suggests that the center-annulus travel 
time difference [6t i] in the separation range A = 10-24° is insensitive to the 
horizontally diverging flow from the centers of the supergranules and should lead 
to a constant signal from the vertical flow. Our measurement of this quantity, 
5.1 ± 0.1 seconds, is constant over the distance range. This magnitude of signal 
cannot be caused by the level of upfiow at cell centers seen at the photosphere of 
10 ms" 1 extended in depth. It requires the vertical flow to increase with depth. 
A simple Gaussian model of the increase with depth implies a peak upward flow 
of 240 ins -1 at a depth of 2.3 Mm and a peak horizontal flow of 700 ms -1 at a 
depth of 1.6 Mm. 
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1. Introduction 



Supcrgranulation, first studied more than 50 years ago QHart, 19 54 Lcig hton, Noyes, a nd Simon, 1962), 



continues to be an active area of research. A comprehensive review ( |Rieutord and Rincon, 2010[ )[ 
details the recent progress. It has proven difficult to measure the vertical flow of 
supergranules at the photospheric level. The recent measurements of a 10 ms -1 
upflow at the center of the average cell with a horizontal variation consis- 
tent with a simple convection cell flDuvall and Birch, 2010[ ) may have finally 
settled the issue. It is only through helioseismology that we would hope to mea- 
sure the supergranulation flows below the photosphere. IDuvall et al.\ (|1997[) first 
showed that time-distance helioseismic techniques are sensitive to supergranules 
and that inversions to derive three-dimensional flow fields might be derived. 
However, I Zhao and Kosovichevl (|2003p showed that their ray-theory inversions 
could not separate the horizontal and vertical flows for a model flow field. 
The use of radial-order filters and Born approximation kernels has led to more 



successful separation flJackiewicz, Gizon, and Birch, 2008 Svanda et at, 2011 1 



|Gizon, Birch, and Spruft] (|2010p summarize the local helioseismic contributions 
to supergranulation. 

An important quantity in time-distance helioseisomology is the arc separa- 
tion [A] between pairs of photospheric locations whose signals are subsequently 
temporally cross correlated. Previous analyses used A less than 5° (61 Mm) 
([Zhao and Kosovichev,~2003| or 2.4° (29.2 Mm) flJackiewicz, Gizon, and Birch, 2008) 1 



Svanda et al, 2011 ), as the sensitivity is greater for the horizontal flow, the 
signal-to-noise ratio is better for small separations, and basically no super- 
granulation signal could be seen at larger A in the travel time difference maps 
constructed from the relatively short (8-12 hours) time intervals required to 
study the one-day lifetime supergranules. In the present work, we show that the 
large separations of A = 10-24° yield center-annulus travel time differences 
[St Q i = outward-going time minus inward going time] from the centers of average 
supergranules that are insensitive to the diverging horizontal flow and hence yield 
a purely vertical flow response. 

In Section [2] the modeling efforts, including flow models that satisfy mass 
conservation, travel times computed from ray theory, and linear wave simulations 
through flow fields are presented. In Section [21 the analysis of data is presented, 
and in Section @] the results are discussed. 



2. Modeling 
2.1. Flow Model 

As our observational technique is to make averages about centers of super- 
granules, the flow model that we take is azimuthally symmetric and decays 
exponentially away from the cell center to simulate the effect of averaging a 
stochastic flow field. A time-independent flow is assumed with axisymmctry 
about the location of the average cell. A Cartesian coordinate system with z 
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positive upwards is denned. For these assumptions, the vanishing divergence of 
the mass flux implies 

d 

Q^(P V *) = -P^h ■ v h , (1) 

where p is the z-dependent density, v z is the z-component of velocity, Vh is 
the horizontal divergence, and Vh is the horizontal velocity. It is assumed that 
horizontal density variations can be ignored, as travel times are only weakly 



sensitive to density changes in the underlying medium (Hanasoge et al, 20121 



The model is assumed to separate into horizontal and vertical functions 

v h = -f(z) S (x,y), (2) 



v z = u(z)V h • g, (3) 

where g(x,y) is a vector in the horizontal plane with no units, /(z) has units 
of velocity, and u{z) has units of velocity times length. Our method of solution 
is to first specify the horizontal function g(x, y) and calculate its horizontal 
divergence Vh • g- The vertical function u{z) is then specified, and Equation (p} 
is used to derive f(z). Some straightforward algebra yields 

n / \ du din/? 

fi*) = Tz +u W-aT- (4) 

In general, models are considered with the horizontal function g(x, y) defined by 

g (x,y)=rJ 1 (kr)e- r / R , (5) 

where f is the outward radial unit vector in a cylindrical coordinate system, J\ is 
the Bessel function of order one, fc is a wavenumber, r is the horizontal distance 
from the origin, and R is a decay length, k and R are, in principle, free parameters 
that will only have the default values k = |^ rad Mm -1 and R = 15 Mm in this 
article. This type of horizontal variation was used by IBirch and Gizonl (|2007p . 
This type of horizontal variation describes well the type of averaging of the 
supcrgranular field done in this paper, where one sees the outward flow as the 
first positive lobe of the J\ function and on average, the inflow from the adjacent 
supergranular cells as the first negative range of the J\ function. 

In general, models have been considered with a Gaussian z-dependence. v z is 
specified at r = as a simple Gaussian: 

v z {r = 0) = ku{z) = v e~ (z - zo)2/2c7 ', (6) 

where zq is the location of the peak of the vertical flow, a z is the Gaussian sigma, 
and vq is the maximum vertical flow. To ultimately explain the large-distance 
travel times, the photosphere needs to be in the far tail of the Gaussian. As 
the upward flow at cell center is 10 ms -1 at the photosphere, this implies a 
considerably larger vertical flow at depth for the average cell. Values for these 
parameters that approximate the data are vq = 240 ms -1 , zq = —2.3 Mm, and 
a z = 0.912 Mm. 
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Figure 1. Comparison of various ray models. The blue curves are for the Gaussian model 
with zo = —2.3 Mm. The lower dot-dashed blue curve is for the horizontal flow component, the 
upper dot-dashed blue curve is for the vertical flow component and the solid blue curve is for 
the sum (upper solid curve). The upper dashed (green) curve is for the sum for zo = —1.15 Mm 
model and the lower dashed (red) curve is for the sum of the zq = —3.45 Mm model. The (lower) 
solid black curve is the sum for a model with a constant upflow of 10ms -1 . 



2.2. Ray Calculations 

In the ray theory the travel time difference for the two directions of propagation 
through a flow v is 

/ r ^, (7) 

where T is the ray path, c is the sound speed, and ds is the element of length in 
the direction of the ray ( |Giles, 2000[ ) . This equation naturally separates into two 
terms for the horizontal and vertical flow contributions to the travel time. The 
ray generation and raw ray integrations are performed using code developed 
and discussed in detail by IGilesI (|2000|) . This code was extended to integrate 
quadrant and annulus surface-focusing integrations for 3D flow and sound speed 
perturbation models. In Figure[T]the black curve shows the travel time difference 
[St i] for a flow model with a constant upflow of 10ms" 1 . The black curve is 
for the sum of horizontal and vertical flows. That is, a horizontal flow of the 
form of Equation ((5|) consistent through the continuity equation with a 10 ms -1 
upflow is used. This signal, with a maximum of 0.6 seconds, is much too small 
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Table 1. Model details. The three models gl, g2, g3 are ones with Gaussian (g) 
upward flow at cell center peaking at zOMm with peak flow v z (zO) ms~ 1 and with 
Gaussian width cr z Mm. The upward flow at the top (z=0) is v z (0)ms~ 1 . The 
peak horizontal flow at the top is ti™ ax (0)ms which occurs at radial distance 
r max Mm. The peak outward flow «™ ax ms _1 occurs at 2(u™ ax )Mm. All of the 
maxima of Vh occur at the same r max by construction. The model siml is the one 
for the small features with just a z flow that are used for the flow perturbation 
described in Section 12.31 



Name 


zo 




vz(0) 


» 4 mal (0) 


j,max 


v z (zO) 


..max 
v h 


2(>™ aX ) 


gl 


-3.45 


1.39 


10 


122 


7 


218 


460 


-2.44 


g2 


-2.3 


0.912 


10 


138 


7 


240 


697 


-1.62 


g3 


-1.15 


0.433 


10 


195 


7 


340 


1609 


-0.83 


siml 


-2.3 


0.816 


6.3 






338 







to explain the results derived in Section 13.31 This implies that the vertical flow 
must increase with depth from the photospheric value of 10 mr 1 . 

Also in Figure [TJ the separate horizontal and vertical travel time contribu- 
tions and the sum are shown for the Gaussian model mentioned above. For this 
relatively shallow model, it is seen that in the distance range 10-24°, the signal 
is essentially constant and is mostly due to the vertical flow. Towards shorter 
distances, the horizontal flow contributes an increasing amount. The Gaussian 
model has three free parameters: vq, zq, and a z . In a later section the travel time 
difference for 10-24° will be determined to be 5.1s. This travel time difference 
and the 10 ms -1 upflow at the photosphere serve to determine two of the three 
free parameters. The travel time at A = 24.5° is forced to 5.1 seconds. Letting 
zq be the free parameter, three models are plotted in Figure [JJ with values of 
z = -1.15 Mm, z = -2.30 Mm, and z a = -3.45 Mm. The shallowest model 
shows the least contribution from the horizontal component in the distance range 
5- 10°. This distance range may supply a way to distinguish the best model while 
comparing with data. 

A summary of the characteristics of these models is contained in Table [TJ The 
maximum vertical flow is at least a factor of 20 larger than the photospheric 
value. The ratio of the maximum horizontal flow to the photospheric value varies 
between four and eight. The maxima of the horizontal flow occurs somewhat 
shallower than the peak vertical flow. 

2.3. Simulation Results 

The applicability of ray theory to time-distance helioseismology has been called 

into question ( |Bogdan, 1997D . There are cases for supergranulation studies (jBirch and Gizonl | 

120071 Figure 5) where it works reasonably well (amplitude of ray theory 25% 

higher than Born approximation kernels) and others where it is discrepant 

by a factor of two (jBirch and Gizonl 120071 Figure 6). One case that works 

extremely well is the comparison of interior rotation determined from global 

modes and time-distance inversions ([Giles , 2000. Figure 6.2) in the radius range 

0.89 < r/ii© < 0.999. The difference between these cases seems to be that the ray 



SOLA: paper.tex; 26 July 2012; 0:29; p. 5 



T.L.Duvall Jr.,S.M.Hanasoge 



8 



S 4 




A [deg] 

Figure 2. Comparison of the center-annulus travel time differences [<5t ;] from the linear 
simulation (blue) with the travel-time difference computed from the ray theory with the 
same flow perturbations (red). The error bars are computed from the scatter far from the 
feature locations. No filtering has been done before the travel time measurements. The average 
travel-time difference in the range A = 10-24° has been scaled to match the observationally 
determined mean 5.1 seconds. The same scaling factor is then used to scale the ray theory 
results. In the range A = 10 — 24°, the ray theory predicts a travel— time difference too large 
by 24%. 



theory works better when the spatial variation of the perturbation is larger, that 
is, much larger than the acoustic wavelength. One way to check the ray theory 
for a particular case is to compare with travel times computed from Born kernels 
(Bir ch and Gizon, 2007[ ). Another way to check is to perform simulations on a 
convectively-stabilized solar model with specified flow velocity perturbations and 
to propagate acoustic waves through the model ( |Hanasoge, Duvall, and Couvidat, 2007| ) .| 
Travel times are computed from the simulations and then compared with ray- 
theory calculations of the travel times using the same flow velocity perturbations. 
This is the type of checking of the ray theory adopted here. 

A global simulation of wave propagation ( Hanasoge and Duvall, 2007 1 with 
768 x 384 grid points in longitude and latitude respectively is performed over 
a ten-hour period. The flow perturbation consists of Gaussian features with 
a constant direction of vertical flow centered at a depth of zq = —2.30 Mm 
with a z = 0.82 Mm and with horizontal er^ = 5.1 Mm. There are 500 of these 
features placed at random longitude-latitude pairs at that depth [zq]. Center- 
annulus travel time differences arc measured from the simulation and averaged 
about the known locations. Realization noise is removed to first order by doing a 
similar simulation with no flow perturbation but with the same source excitations 
QWerne, Birch, and Julien,~2 004; Hanaso ge, Duvall, and Couvidat, 2007| ) and sub-| 
tracting the resultant 5t j to obtain noise-corrected results. 

The flow model that is inserted into the simulation is used with Equation ([7]) 
to derive ray-theory predictions of the travel time differences. The results of 
the ray-theory computations and the travel time difference measurements from 
the simulation are shown in Figure [5] As the simulation only contains vertical 
flows, it is expected that there would be little or no variation of the travel-time 
differences with A. The ray theory does predict too large a travel time difference 
by 24%. This excess amplitude is similar to that found by IBirch and Gizonl 
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(2007). Although the amplitude computed from the ray theory is not precise, the 
necessity of large vertical flows to generate a travel-time difference of 5.1 seconds 
at A = 10-24° is confirmed. The peak vertical flow for the normalized case is 
338 ms -1 . The model parameters are detailed in Table [TJ 



3. Data Reduction and Analysis 

3.1. Reduction 

Dopplcrgrams from the Helioseismic and Magnetic Imager (HMI) : ( |Schou et al., 2012[ )[ 
onboard the Solar Dynamics Observatiory (SDO) spacecraft were analyzed for 
the present work. 32 days (10 June- 11 July 2010) of Dopplcrgrams were used 
to derive the final results, although a number of tests subsequently mentioned 
used only the last three days, 9-11 July 2010. This particular time period was 
used as the Sun was very quiet (sunspot number RI w 15 for these two months), 
and it was the final two-month continuous coverage period for the Michelson 
Doppler Imager (MDI) instrument (Scherre r et al., 1995| ). It might be useful to 



compare the results from the two instruments, but only the HMI data was used 
for the present study. 

Raw Dopplergrams have a nine-day averaged Dopplergram subtracted as well 
as the spacecraft velocity which has a significant 24-hour component. These 
corrected Dopplergrams are remapped onto a coordinate system with equal 
spacing in latitude and longitude of 0.03° covering a range of 144° in both 
latitude and longitude. The remapping is achieved using a bilinear interpolation. 
The remapping in longitude is onto the Carrington system at a central longitude 
that crosses central meridian at the middle of the twelve-hour interval used. Two 
twelve-hour intervals are used for each day, covering the first and second halves 
of the day. Each remapped image is Fourier-filtered and resampled at 0.24° 
per pixel. Twelve-hour datacubes are constructed from these individual images 
with 600 x 600 spatial points and 960 temporal points for the 45 s sampling. 
This procedure works well for the large A emphasized in the present study, but 
would need to be modified to study smaller A. Studying smaller A is deferred 
until future work. 

The centers of the supergranules are located by a procedure used previously 



(Du vall and Birch, 2010 ). The datacubes are spatio-temporally filtered to pass 
just the solar /-mode oscillations. Center-annulus travel time differences are 
computed ( |Duvall and Gizon, 2000[ ) for the distance range 0.48-1.02°. By us- 
ing the difference of inward-going times minus outward-going times, the maps 
are equivalent to maps of the horizontal divergence. The equivalence of these 
maps to maps of supergranules has been shown before ( |Duvall and Gizon7"2 000; 



Duvall and Birch, 2010 1. The travel-time difference maps arc smoothed with a 
Gaussian of a = 2.9 Mm to more easily determine the cell-center locations. Local 
maxima of this signal are picked as the cell centers. Lists of the locations are 
made. The one of a pair with the smaller divergence that are closer together 
23 Mm is rejected. As it is desired to use center-annulus separations A <= 25°, 
locations close to the edge of the maps are not used. On average, 930 features 
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Table 2. Nominal phase speed filter parameters. The first col- 
umn is an identification number. The second column is the range 
in distance [A]. The third column is the mean distance [A mea n]- 
The fourth column is the central phase speed. The fifth column is 
the full width at half maximum (FWHM) of the Gaussian phase 
speed filter. For the tests in Section 13.21 filters 1-13 were used 
while for the final results in Section 13.31 filters 1-14 were used. 



Filter 


A[°] 


A mean [ ] 


^0 [ km s 1 ] 


<5u[kms - 1 ] 


1 


1.56-2.28 


1.92 


32.0 


13.2 


2 


2.28-3.00 


2.64 


41.8 


4.2 


3 


3.00-3.96 


3.48 


46.9 


6.0 


■1 


3.96-4.92 


4.44 


54.7 


9.5 


5 


4.92-6.12 


5.52 


64.5 


10.1 


6 


6.12-7.80 


6.96 


75.1 


11.1 


7 


7.56-9.24 


8.40 


84.3 


10.4 


8 


9.12-10.7 


9.84 


93.5 


9.4 


9 


10.92-12.60 


11.76 


104.4 


9.7 


10 


12.36-14.52 


13.44 


114.2 


12.5 


11 


14.28-16.44 


15.36 


125.2 


12.3 


12 


17.16-19.32 


18.24 


141.9 


12.6 


13 


19.08-22.20 


20.64 


156.0 


18.5 


14 


22.08-24.00 


23.04 


170.4 


11.8 



are located for each datacube with a total of 59549 for the 64 twelve-hour 
datasets. Overlays of cell-center locations with the divergence maps are shown 



in (Duvall and Birch, 2010). 



3.2. Filtering 

Ray theory has proven invaluable in the development of time-distance helioseis- 
mology ( |Duvall et ai, 1993] ) . It led to the basic idea that signals propagating 
along ray paths between surface locations would lead to correlations between 
the temporal signals from which travel times could be inferred. As acoustic, or 
p, modes with the same horizontal phase speed \uj/k: tu angular frequency and 
k horizontal wavenumber] travel along the same ray path to first order, it was 
natural to filter the data in ui/k to isolate waves traveling to a certain depth. This 
type of filter was introduced by IDuvall et al.\ jTT>').''»! with subsequent develop- 
ment bv lKosovichev and Duvalll(|1997p and [Giles ( 20DS1)- I n the spatio-temporal 
power spectrum of solar oscillations {k~uj diagram), a line from the origin has a 
constant phase speed uj/k. A range of phase speeds thus corresponds to a pie- 
shaped wedge. Phase-speed filtering then corrcponds to multiplying the Fourier 
transform of the input datacube by a pie-shaped wedge (appropriately modified 



for 3D) with some Gaussian tapering in phase speed ( |Couvidat and Birch, 2006 1 . 

The overriding principle of the use of the phase-speed filter is that solar per- 
turbations lead to refractive time changes for signals travelling along a ray path. 
That this model is deficient was shown nicely in the work offCouvidat and Birchl 
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Figure 3. Center-annulus travel time differences [<54 ;] averaged about the supergranulc 
centers for different filters for the thirteen A ranges for three days of HM1 data. From the 
bottom to top (on the right side), the black (x) points are for the the nominal phase-speed 
filters. The blue (*) points are for the phase-speed filters with the widths doubled. The red (o) 
points arc for the width tripled. The orange (square) points are for a filter width that is half the 
nominal phase speed. The magenta (diamonds) curve is for the constant dcgrec-width-filtcr 
(Section 13. 311 with width Vi = 400. The green (pentagrams) curve is for no phase-speed filter. 
In all cases the / mode is excluded via filtering as is signal outside the frequency bandpass 
1.5 < v < 6mHz. These results were obtained using the three days of data 9-11 July 2010. 

(2006) in which an artificial signal was detected with different amplitude de- 
pending on the width [5v] of the pie-shaped wedge. They concluded that the 
Born-approximation kernels ( |Birch and Kosovichev, 2000[ ), which are derived 
based on a single scattering from solar perturbations, arc more appropriate. 
A better way to think of the problem is to consider waves impinging on a 
perturbation which are subsequently scattered. 

An initial analysis was performed using some nominal phase speed filters and 
distance ranges detailed in Table [2J Thirteen distance ranges (and hence thirteen 
different filters) were applied to the six input datacubes for the dates 9-11 July 
2010. The overall distance range covered for this test is 3.0-22.2°. Center- 
annulus travel-time differences [<5i i] were obtained for the thirteen distance 
ranges for each of the six datacubes using a standard surface-focusing time- 
distance analysis with the Gizon-Birch method of extracting travel times from 
the correlations ( |Gizon and Birch, 2004[ ). The travel times are averaged over 
the distance range and about the locations of the supcrgranular centers and the 
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signal at the average cell center location is extracted. Uncertainties are estimated 
from the scatter of the six values. The travel-time differences are plotted versus 
A in Figure [3] (black). The variation with A does not agree with the models in 
Figurc[TJ but peaks near A = 5° and decays approximately to zero near A = 20°. 

There was some concern that this distance dependence might have something 
to do with the filtering, so a case was analyzed with no phase-speed filtering. 
The / mode was excluded and there was a frequency bandpass filter transmitting 
v = 1.5-6mHz. The results of that analysis are also shown in Figure [3] (green). 
This shows a behavior with A much closer to the models in Figure [1] but with 
considerably more noise and a much larger signal close to five seconds. It was 
concluded that the phase speed filters were significantly degrading the signal. 
To test this further, the filters were broadened by the factors 2, 3, 4, 5 with the 
results of factors 2 (blue) and 3 (red) also shown in Figure [3] More signal is 
obtained with the larger-width filters, but even with the factor of 5, only about 
50% of the signal is obtained and there is still a decay of the signal towards 
larger A. 

Because of the approximate constancy of the simulation travel times [St a i] 
with A in Section 12.31 and the approximate constancy with A for the no phase- 
speed filter case, it was concluded that the decay of 5t Q i toward larger A in 
the phase-speed filtered results is an artifact of the filtering. Filtering would 
seem to be desirable because of the large difference in errors for the unfiltered 
and filtered cases. But a filter is needed that gives a more unbiased result at the 
large A in order that the results not be overly dependent on the modeling. A first 
attempt at this was to use a set of filters with the full width at half maximum 
(FWHM) equal to half the phase speed ( |Birch et ai, 2006 ) . The results of these 
measurements are also shown in Figure [3] (orange). The dependence on A is 
somewhat reduced but there is still some decay at the largest A. 

A different way to look at this problem is to consider a p mode with frequency 
v and spherical harmonic degree £ impinging on a supcrgranulc. To the extent 
that the supcrgranulc is static, the mode frequency v will be maintained and 
the scattering will spread power in the spectrum away from the nominal value 
of I. The scattering should spread power over a width in I that depends on 
the spectrum of supcrgranulation, which peaks near I = 120. For similar values 
of the phase speed (which can also be characterized by v/i), the spread power 
should be over a similar width in i QChou et al., 1996} |Woodard, 2002[ ). A filter 
that is centered at a certain phase speed but whose width is constant with v 
and characterized by the FWHM in I (= Ti) may be what is needed. A value 
of Yi — 400 should capture most of the supergranular signal. Such a filter was 
implemented and applied to the thirteen distance ranges used. The results are 
also shown in Figure [3J The approximate constancy with A > 5° and the larger 
values than for the normal phase speed filters suggests that this type of filtering 
is a useful concept. The filter used has a flat top of width T^/2 and is tapered 
to zero by a cosine bell that goes from 1 to in Tg/2. 

Some additional tests of this new filter concept were computed. In Figure[UJa) 
is shown the peak supergranular signal for A = 19.08-22.2° as a function of the 
filter width Ti. For small Ti, much of the supergranular signal is removed, but 
as the width is increased more and more of the signal is transmitted. There 
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Figure 4. (a) Travel-time difference [<5t i] versus filter FWHM [I^]. The unfiltered case has 
5.3 ± 1.2 seconds, (b) The travel time difference from (a) divided by the size of the error bar 
from (a) versus the filter FWHM Tp. The value for the unfiltered case is 4.6. 



is clearly a variation in the size of the error bars across Figure IHa). Another 
thing to examine is the ratio of the signal in Figure EJa) to the size of the error 
bar, which should give a value of the signal to noise ratio (S/N) of the filter 
( Couvidat and Birch, 2006 ). This is plotted in Figure E£b). A rather broad peak 
is seen more or less centered on = 400. One noticeable aspect of Figure HJb) 
is that the unfiltered case gives a worse S/N (4.6) than the filters near the peak. 
Also narrow filters have reduced S/N. Based on the flatness of the Ti = 400 
results with A in Figure [3J the extraction of essentially all the supergranular 
signal in Figure Eta), and the value of the S/N in Figure E^b), Ti — 400 has 
been chosen for further analysis. 



3.3. T e = 400 Analysis 



The 64 twelve-hour datacubes were analyzed with the fourteen phase-speed 
filters with = 400. Center-annulus cross correlations were computed for 
the various distance ranges of Table 1. For checking purposes, it would be 
useful to be able to compute travel times using both the Gabor wavelet method 
dDuvall et al, 1997| ) and the Gizon-Birch method ( jGizon and Birch, 2004| |. How-| 
ever, for the large distances in the present study (A < 24°), individual twelve- 
hour cross correlations have insufficient signal-to-noise ratio to enable the Gabor 
wavelet fitting. It was decided to average the temporal correlations spatially 
about the supergranular centers. The averaging of the 930 (on average) corre- 
lations yields a high signal-to-noise ratio correlation that is amenable to the 
Gabor wavelet fitting at the largest distances [A]. Averaging the resultant <5i i 
for the 64 twelve-hour datacubes enables the use of the envelope-time differences 
[tenv] as well as the normal phase times [t p h]- Although most theoretical work 
applies to phase time differences, it will be useful to examine the envelope-time 
differences as they may have independent information. 

It was found previously that there is an offset from zero (of unknown origin) 
of the center-annulus travel-time differences ( |Duvall and Gizon, 2000[ ) . For the 
small A of that study, the value of the offset was 0.16 ± 0.02 seconds. This offset 
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Figure 5. (a) Center-annulus travel-time difference averaged about supergranule centers for 
the range of A = 22.08-24°. The scale of the colorbar at right is in seconds. Only the central 
point corresponding to the supergranular center is used in the present study, (b) Azimuthal 
average of (a). Note the offset at large radii. This offset is believed to be an artifact which 
needs to be removed from the results, (c) The offset at r = 58 Mm for the different travel-time 
definitions versus A. Blue is for the Gabor wavelet phase time differences. Red is for the 
Gizon-Birch phase time differences and the green is for the Gabor wavelet envelope time 
differences, (d) The resultant travel time differences averaged for the 64 12-hour datacubcs 
corrected for the offset in (c). The colors are the same as in (c). 



is assumed unrelated to the flows to be measured and needs to be removed 
from computed travel times. To measure it for the present study, the average 
maps about the supergranule centers for the 64 12-hour datacubes are com- 
puted for the different distance ranges. An example is shown for the largest 
range of A = 22.08-24° in Figure [UJa). An azimuthal average of Figure [5ja) 
is computed and the result is shown in Figure [SJb). At large radii [r], the the 
signal is approximately constant. This constant, which is assumed to be the 
offset that needs to be subtracted from the results, is plotted for the various A 
and travel-time methods in Figure 03c) . The Gizon-Birch and Gabor wavelet 
phase times yield almost exactly the same results, which is expected for the 
quiet Sun and the use of the same time window for fitting (20 minutes). The 
travel-time differences corrected for this offset are plotted in Figure E^d). Again, 
the Gizon-Birch and Gabor wavelet phase times yield almost precisely the same 
results with an average of 5.1 ±0.1 seconds on the range A = 10 - 24°. The Gabor 
wavelet envelope times yield a considerbly larger value of 9.7 ± 0.3 seconds on 
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Figure 6. Comparison of the three ray models from Figure^ with the HMI results (black 
with symbols and errors) of the = 400 filtering and the Gabor-wavelet phase— speed time 
differences from Figure (3d). These models are the sum of vertical and horizontal signals for 
the three Gaussian vertical flows peaking at zo = —3.45 Mm (red dot— dashed), zo = —2.3 Mm 
(blue solid), and zq = —1.15 Mm (green dashed) specified in Table [T] 



A = 10-24°. In general the group velocity, duj/dk is about one half of the 
phase speed u>/k for p modes, possibly leading to the larger travel time. The 
systematic error from Figure [5jc) is generally less the 10% of the signal shown 
in Figure (5Jd). 

The Gabor-wavelet phase-time differences are plotted in Figure [S] versus the 
three ray models of Figured] The increase of the observed times for A = 2.5-8° 
may be the effect of the horizontal flow. The agreement is a little better for the 
zo = —2.3 Mm model. 



4. Discussion 

The main observational result of this article is that in the distance range A = 10 - 
24° that the mean travel-time difference [5t Q i] at the center of the average 
supergranule is 5.1 ± 0.1 seconds. How secure is this result? There are no other 
results to compare with in this distance range. Also, the method of averaging over 
many supergranule s is uncommon ( |Birch et ai, 2006} |Duvall and Birch, 2010| 



Svanda et aL, 2011). The largest distance [A] used previously is 5° bv lZhao and Kosovichev[ 
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(|2003p . The phase-speed filter used in that study would have reduced the travel- 
time difference [<5t i] by a factor of five and so it is difficult to compare. In 
addition, only inversions were presented and not raw travel-time differences. 

One of the most interesting aspects of the present results is the large factor 
(> 20) between the photosphcric vertical flow and the peak vertical flow. This 
increase of vertical flow then also requires an increase of the horizontal flow from 
the photosphere to the peak of a factor between 3.7 and 8.3 from the simple 
Gaussian models gl-g3. Some simple tests suggest that this horizontal velocity 
should be detectible by a quadrant time-distance analysis, which will be left for 
future work. Some previous analyses have detected some larger flows below the 
surface than at the surface ( |Duvall et ai, 1997 ). Possibly the subsurface flow 
is highly variable which has made average properties difficult to ascertain from 
inversions of individual realizations. 

These results confirm that supergranulation is a shallow phenomenon, with 
the flow peaking within a few Mm of the photosphere. This would tend to 
support models in which supergranulation is related to a near-surface phe- 
nomenon such as granulation (Rast, 2003). The supergranular wave measure- 



ments ( Gizon, Duvall, and Schou, 2003 Sc 



rou, 2003 ), and in particular the vari- 



ation with latitude of the anisotropy of wave amplitude, would suggest a con- 
nection with rotation. The near-surface shear layer has been modeled as a 
possible exciter of supergranulation ( Green and Kosovichev, 2006[ ). |Hathaway| 



(1982) found that supergranulation as a convective phenomenon generates a 
near-surface shear layer. A difficulty with the near-surface shear layer as su- 
pergranular exciter is that recent simulations of IStein et al.\ (|2009[) develop a 
near-surface shear layer but no hint of excess power at supergranular scales is 
seen. 
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